import matplotlib.pyplot as plt
import numpy as np
from sklearn.datasets import load_iris

# 加载鸢尾花数据集
x, y = load_iris(return_X_y=True)

# 1. 数据中心化：减去每个特征的均值
x_centered = x - np.mean(x, axis=0)

# 2. 计算协方差矩阵
cov_matrix = np.cov(x_centered.T)  # 转置以使特征为行

# 3. 特征分解：计算协方差矩阵的特征值和特征向量
eigenvalues, eigenvectors = np.linalg.eig(cov_matrix)

# 4. 选择主成分：按特征值降序排序，取前两个特征向量
sorted_indices = np.argsort(eigenvalues)[::-1]  # 降序排列索引
eigenvectors_sorted = eigenvectors[:, sorted_indices]
components = eigenvectors_sorted[:, :2]  # 取前两个主成分

# 5. 数据投影：将中心化数据投影到主成分上
reduced_x = np.dot(x_centered, components)

# 按类别分类降维后的数据点
red_x, red_y = [], []    # 存储类别0（setosa）的x、y坐标
blue_x, blue_y = [], []  # 存储类别1（versicolor）的x、y坐标
green_x, green_y = [], [] # 存储类别2（virginica）的x、y坐标

for i in range(len(reduced_x)):
    if y[i] == 0:
        red_x.append(reduced_x[i][0])
        red_y.append(reduced_x[i][1])
    elif y[i] == 1:
        blue_x.append(reduced_x[i][0])
        blue_y.append(reduced_x[i][1])
    else:
        green_x.append(reduced_x[i][0])
        green_y.append(reduced_x[i][1])

# 绘制散点图可视化结果
plt.scatter(red_x, red_y, c='r', marker='x', label='Setosa')
plt.scatter(blue_x, blue_y, c='b', marker='D', label='Versicolor')
plt.scatter(green_x, green_y, c='g', marker='.', label='Virginica')
plt.xlabel('First Principal Component')
plt.ylabel('Second Principal Component')
plt.legend()
plt.title('PCA Projection of Iris Dataset (Manual Implementation)')
plt.show()